/*******************************************************************************
	- Generate figures in the main manuscript
*******************************************************************************/


********************************************************************************
* initial settings
********************************************************************************

version 14.2
clear 
drop _all
set more off

*----------------------------------

// set your own local path where you put the JAERE_Submission folder
if c(username) == "zhenx" {
	global LOCALPATH "C:/Users/zhenx/Dropbox/Research/Temperature and Athletes/temp_acclimatization/Code" 
}
global PATH "$LOCALPATH/JAERE_Submission"
global REG "$PATH/Results/Regressions/Baseline"
global FIGURE "$PATH/Results/Figures/Baseline"
global DATA "$PATH/Data"

cd "$FIGURE"

********************************************************************************
* Figure 1: nonlinear relation between temperature and performance
********************************************************************************

*****************************coefficient estimates******************************

import delimited "$REG/Temperature_Effect_by_Event_Type.txt", clear

input 
	60 0 0 0
end

gen j = "coef" 	 if v1 != ""
replace j = "se" if v1 == ""
carryforward v1, replace

drop in 1/1
drop if strmatch(v1,"*Constant*") | strmatch(v1,"*Obser*") | strmatch(v1,"*R-squa*") | ///
		strmatch(v1,"*Adjusted*") | strmatch(v1,"*Numbe*") | strmatch(v1,"*Robust*") | ///
		strmatch(v1,"*p<0.01*")
drop if v2 == ""
nrow

ren coef j
ren VARIABLES var
replace var = subinstr(var,".Meet_temp_bin","",.)
destring var, replace

label define var_label 40 "<40" 45 "40-45" 50 "45-50" 55 "50-55" 60 "55-60" 65 "60-65" 70 "65-70" 75 "70-75" 80 ">75"
label values var var_label

global varlist Sprint Strength Endurance
foreach v in $varlist {
	replace `v' = subinstr(`v',"*","",.)
	replace `v' = subinstr(`v',"(","",.)
	replace `v' = subinstr(`v',")","",.)
	destring `v', replace
}

reshape wide $varlist, i(var) j(j) string
foreach v in $varlist {
	gen `v'H = `v'coef + 1.96 * `v'se
	gen `v'L = `v'coef - 1.96 * `v'se
}

grstyle init
grstyle set plain, nogrid compact

foreach v in $varlist {
	twoway ///
		(connected `v'coef var, lcolor(navy) lwidth(medthick) msize(medium) mcolor(navy) msymbol(circle)) ///
		(line `v'H var, lcolor(gs10) lwidth(medthick) lpattern(dash) cmissing(n)) ///
		(line `v'L var, lcolor(gs10) lwidth(medthick) lpattern(dash) cmissing(n)) ///
		, ///
		xlabel(40(5)80, valuelabel) ///
		ylabel(-0.2(0.05)0.05) ///
		yline(0, lpattern(solid) lcolor(red)) ///
		legend(size(small) order(1 "Estimate" 2 "95% CI") position(6) ring(0) region(fcolor(none) lcolor(none))) ///
		title(`v') ///
		xtitle("") ytitle("Performance Change Relative to 55-60F") ///
		xsize(6) ysize(6) 
	graph save "Temp_`v'.gph", replace
}

***************************temperature distribution***************************

*** import data 
	use "$DATA/TF_Data.dta", clear

*** drop training temperature lower than 40F to mitigate noise
	/*avoid the case where athletes might get trained indoors during extreme cold weathers*/
	drop if Team_tavg <= 40

*** drop the observations that have missing values on key variables	
	drop if Meet_tavg == . | Team_tavg_7 == . | Meet_preciptotal ==. | Meet_ozone == .
	
*** keep variables of interest 
	keep event_type Meet_tavg

*** plotting
	foreach event in Sprint Strength Endurance {
		histogram Meet_tavg if Meet_tavg >= 40 & Meet_tavg <= 75 & event_type == "`event'" ///
			, bin(45) fcolor(gs12) lcolor(gs14) xlabel(40(5)75) xsize(6) fysize(25) xtitle(Mean Daily Temperature (F)) ///
			saving("Meet_Tavg_Dist_`event'.gph", replace)
	}

********************************combined graph********************************

foreach v in Sprint Strength Endurance {
	graph combine ///
		"Temp_`v'.gph" "Meet_Tavg_Dist_`v'.gph" ///
		, col(1) imargin(0 0.6 0 0) xsize(10) scale(1.2) ///
		caption(a, size(huge) position(10) margin(zero))
	graph save "Temperature_Effect_by_Event_`v'.gph", replace
}

graph combine ///
	"Temperature_Effect_by_Event_Endurance.gph" ///
	"Temperature_Effect_by_Event_Sprint.gph" ///
	"Temperature_Effect_by_Event_Strength.gph" ///
	, col(3) imargin(0 0.6 0 0) xsize(10) scale(1.2)
	gr_edit .plotregion1.graph2.caption.text = {}
	gr_edit .plotregion1.graph2.caption.text.Arrpush b
	gr_edit .plotregion1.graph3.caption.text = {}
	gr_edit .plotregion1.graph3.caption.text.Arrpush c
graph export "Temperature_Effect_by_Event.eps", replace

foreach v in Sprint Strength Endurance {
	erase "Temp_`v'.gph"
	erase "Meet_Tavg_Dist_`v'.gph"
	erase "Temperature_Effect_by_Event_`v'.gph"
}

********************************************************************************
* Figure 2: nonlinear relation between temperature and performance 
********************************************************************************

*** generate a dataset to plot the figure using estimation results
	clear
	set obs 1000
	gen n = _n
	gen Meet_tavg = 30 + n*(90-30)/500 if n <= 500
	replace Meet_tavg = 30 + (n-500)*(90-30)/500 if n > 500
	gen hot = cond(n<=500,0,1)
	keep if Meet_tavg >=40 & Meet_tavg <= 80

*** set other control variables at the mean level
	gen wind = -0.0000445
	gen Meet_preciptotal = 0.0819
    gen Meet_ozone = 0.0367131
    gen travel_dist = 277.2341

*** create splines
    mkspline2 meet_cubspl= Meet_tavg, cubic di knots(49.6 56.7 63.8)
	mkspline meet_linspl1 50 meet_linspl2 60 meet_linspl3  = Meet_tavg, di

*** plot the figure 
	foreach model in Quadratic Cubic Cubic_Spline Linear_Spline {
		preserve
			estimates use "$REG/Endurance_`model'"
			predict yhat, xb
			predict sehat, stdp
			replace yhat = yhat 
			gen yhat_H = yhat + 2*sehat 
			gen yhat_L = yhat - 2*sehat 

			
			grstyle init
			grstyle set plain, nogrid
			twoway ///
				 (line yhat_H Meet_tavg if hot == 1, lwidth(medium) lcolor(black) lpattern(dash) cmissing(y) sort) ///
				 (line yhat_L Meet_tavg if hot == 1, lwidth(medium) lcolor(black) lpattern(dash) cmissing(y) sort) ///
				 (line yhat   Meet_tavg if hot == 1, lwidth(medium) lcolor(black)  cmissing(y) sort) ///
				 (line yhat_H Meet_tavg if hot == 0, lwidth(medium) lcolor(gs10) lpattern(dash) cmissing(y) sort) ///
				 (line yhat_L Meet_tavg if hot == 0, lwidth(medium) lcolor(gs10) lpattern(dash) cmissing(y) sort) ///
				 (line yhat   Meet_tavg if hot == 0, lwidth(medium) lcolor(gs10)  cmissing(y) sort) ///
				 , ///
				 legend(order(1 "+/- 2SE" 3 "Warmest Pre-Exposure" 4 "+/- 2SE" 6 "Coolest Pre-Exposure") position(6) ring(0) region(fcolor(none) lcolor(none))) ///
				 xtitle("Temperature (F)") ytitle("Standardized Performance Relative to World Record") ///
				 xlabel(40(5)80) ///
				 ylabel(-2.6(0.1)-2.3) /// 
				 title(`model', size(medium)) ///
				 plotregion(margin(zero)) ///
				 scale(1.2) ///
				 caption(a, size(huge) position(10) margin(zero))
			graph save  "Endurance_`model'.gph", replace
		restore
	}

*** combine the figures together
	grstyle set plain, nogrid 
	grc1leg2 "Endurance_Quadratic.gph" "Endurance_Cubic.gph" "Endurance_Linear_Spline.gph" "Endurance_Cubic_Spline.gph" ///
		, legendfrom(Endurance_Quadratic.gph) col(4) xsize(10) ysize(16)
	gr_edit .style.editstyle margin(zero) editcopy
	gr_edit .plotregion1.graph1.style.editstyle margin(medium) editcopy
	gr_edit .plotregion1.graph2.style.editstyle margin(medium) editcopy
	gr_edit .plotregion1.graph3.style.editstyle margin(medium) editcopy
	gr_edit .plotregion1.graph4.style.editstyle margin(medium) editcopy
	gr_edit .legend.style.editstyle margin(vmedsmall) editcopy
	gr_edit .style.editstyle declared_xsize(12.5) editcopy
	gr_edit .plotregion1.graph3.title.text = {}
	gr_edit .plotregion1.graph3.title.text.Arrpush Linear Spline
	gr_edit .plotregion1.graph4.title.text = {}
	gr_edit .plotregion1.graph4.title.text.Arrpush Cubic Spline
	gr_edit .plotregion1.graph2.caption.text = {}
	gr_edit .plotregion1.graph2.caption.text.Arrpush b
	gr_edit .plotregion1.graph3.caption.text = {}
	gr_edit .plotregion1.graph3.caption.text.Arrpush c
	gr_edit .plotregion1.graph4.caption.text = {}
	gr_edit .plotregion1.graph4.caption.text.Arrpush d
	graph export "Adapt_Function_Forms.eps", replace
	
*** delete useless intermediate files
	foreach file in "Endurance_Quadratic.gph" "Endurance_Cubic.gph" "Endurance_Linear_Spline.gph" "Endurance_Cubic_Spline.gph" {
		erase `file'
	}




